#include "exceptions.hpp"
#include "eigenConcat.hpp"
#include "typeHelpers.hpp"
#include <iostream>

namespace hms {


template<typename Derived>
void printAsStructMesh(
	const Eigen::DenseBase<Derived>& vals, Array2i meshSize
){
	static constexpr Index ng {1};
	Index nx, ny;
	if ( vals.cols() == meshSize.prod() ){
		/* domain values */
		nx = meshSize.x();
		ny = meshSize.y();
	} else if ( vals.cols() == (meshSize + 2*ng).prod() ){
		/* all values */
		nx = meshSize.x() + 2*ng;
		ny = meshSize.y() + 2*ng;
	} else {
		std::cout << "printAsStructMesh: Wrong size! "
		<< "\nMesh size: " << meshSize.x() << "x" << meshSize.y()
		<< ", values:\n" << vals
		<< "\n";
		return;
	}
	std::cout << vals.reshaped(nx * vals.rows(), ny)
		.transpose().colwise().reverse() << "\n";
}


template<typename Derived, typename MeshType>
void printAsStructMesh(
	const Eigen::DenseBase<Derived>& vals, const StructMeshBase<MeshType>& mesh
){
	if constexpr (std::is_base_of_v<StructMeshBase<MeshType>, MeshType> )
		printAsStructMesh( vals, mesh.size() );
}


/* this function works in debug mode just fine and produces the correct
 * results. however, in release mode, it causes a segfault. */
template<typename MeshType, typename Derived>
Array2Xi parentCells(
	const StructMeshBase<MeshType>& mesh,
	const Eigen::ArrayBase<Derived>& edgeIndices,
	int order
){
	static constexpr Index ng {1};
	const Index
		nx { mesh.nCellsX() },
		ny { mesh.nCellsY() },
		ni { edgeIndices.size() },
		nh { nHorizontalEdges(nx, ny) },
		nxg { nx + 2*ng };

	static_assert( remove_qualifiers<decltype(edgeIndices)>::ColsAtCompileTime == 1 );

	const auto nxgOr1 { edgeIndices.unaryExpr( [&](Index i){
		return i < nh ? nxg : 1;
	}) };

	const auto rowIndex { edgeIndices.unaryExpr( [&](Index i){
		return i < nh ?
			i / nx :
			(i-nh) / (nx+1);
	})};
	const auto colIndex { edgeIndices.unaryExpr( [&](Index i){
		return i < nh ?
			i % nx :
			(i-nh) % (nx+1);
	})};

	RowArray2i sides {-order, order-1};

	/* with this setup, the cell indices are:
	 * sides * nxgOr1 + nxg * (rowIndex + ng) + colIndex + ng */
	return { (
		sides.replicate(ni, 1) * nxgOr1.template replicate<1,2>()
		+ ( nxg * (rowIndex + ng) + colIndex + ng ).template replicate<1,2>()
	).transpose() };
}

}